;-------------------------------------------------------------LICENSE--------------------------------------------------------------;
;                                                                                                                                  ;
;The MAP code is written in Fortran language for magnetohydrodynamics (MHD) calculation with the adaptive mesh refinement (AMR)    ;
;and Message Passing Interface (MPI) parallelization.                                                                              ;
;                                                                                                                                  ;
;Copyright (C) 2012                                                                                                                ;
;Ronglin Jiang                                                                                                                     ;
;rljiang@ssc.net.cn                                                                                                                ;
;585 Guoshoujing Road. Pudong, Shanghai, P.R.C. 201203                                                                             ;
;                                                                                                                                  ;
;This program is free software; you can redistribute it and/or modify it under the terms of the GNU General Public License         ;
;as published by the Free Software Foundation; either version 2 of the License, or (at your option) any later version.             ;
;                                                                                                                                  ;
;This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of    ;
;MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License for more details.                        ;
;                                                                                                                                  ;
;You should have received a copy of the GNU General Public License along with this program; if not, write to the Free Software     ;
;Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301, USA.                                                   ;
;                                                                                                                                  ;
;-------------------------------------------------------------LICENSE--------------------------------------------------------------;

pro div_b_2d, bx, by, x, y, NG = ng

   if (n_elements(bx) ne n_elements(by)) then print, ' The shapes of bx and by do not conform '
   if (n_elements(bx[*, 0]) ne n_elements(x)) then print, ' The coordinate of x is not correct.'
   if (n_elements(bx[0, *]) ne n_elements(y)) then print, ' The coordinate of y is not correct.'
   if (n_elements(ng) eq 0) then ng = 4

   nx = n_elements(x)
   ny = n_elements(y)

   xx = fltarr(nx, ny)
   yy = fltarr(nx, ny)
   div_b = fltarr(nx, ny)

   for j = 0, ny - 1 do xx[*, j] = x
   for i = 0, nx - 1 do yy[i, *] = y

   div_b[1 : nx - 2, 1 : ny - 2] =                                                                                                 $
      (bx[2 : nx - 1, 1 : ny - 2] - bx[0 : nx - 3, 1 : ny - 2]) / (xx[2 : nx - 1, 1 : ny - 2] - xx[0 : nx - 3, 1 : ny - 2]) +      $
      (by[1 : nx - 2, 2 : ny - 1] - by[1 : nx - 2, 0 : ny - 3]) / (yy[1 : nx - 2, 2 : ny - 1] - yy[1 : nx - 2, 0 : ny - 3])

   shade_surf, div_b[ng : nx - ng - 1, ng : ny - ng - 1], x[ng : nx - ng - 1], y[ng : ny - ng - 1], charsize = 4

end
